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We study a generalized quantum hard-core dimer model on the square and honeycomb lattices, 
allowing for first and second neighbor dimers. At generalized RK points, the exact ground states 
can be constructed, and ground-state correlation functions can be equated to those of interacting 
1 + 1 dimensional Grassmann fields. When the concentration of second neighbor dimers is small, 
the ground state correlations are shown to be short-ranged corresponding to a (gaped) spin liquid 
phase. On a 2-torus, the ground states exhibit fourfold topological degeneracy. On a finite cylinder 
we have found a dramatic even-odd effect depending on the circumference, and propose that this 
can be used as a numerical diagnostic of gapped spin-liquid phases, more generally. 



There has been a surge in numerical experiments [THS] 
on simple models reporting evidence of the existence of 
fully gapped (short-ranged) spin-liquid states, including 
reports for the spin-1/2 Heisenberg model with first and 
second neighbor anti- ferromagnetic couplings (Jj and J2) 
on the square [3] and honeycomb [4] lattice as well as the 
Hubbard model on the honeycomb lattice [2]. For the 
square lattice model, an early numerical study [5] on a 
relatively small system was interpreted as evidence for a 
spin-liquid, but this result was seriously questioned on 
the basis of the results of series expansion [6], large- 7V[7], 
and large- S [8j studies of the same model. In the case 
of the honeycomb lattice, recent studies reported[2] ev- 
idence of such a state for an intermediate range of U/t 
in the Hubbard model. For both the square and honey- 
comb lattice, there is no solvable model with the same 
symmetries as the Heisenberg model for which a gaped 
spin-liquid state has been found. Having a caricature of 
the putative state is useful in thinking about the phase, 
and possibly also in designing numerical experiments to 
confirm or falsify its existence. 

Here, we report a generalization of the quantum hard- 
core dimer model^HO] on the square and honeycomb lat- 
tices, for which the exact spin-liquid groundstates can be 
obtained for certain values of coupling constants dubbed 
as "generalized Rokhsar-Kivelson (RK) points". In the 
limit in which all dimers occupy nearest-neighbor links, 
this spin-liquid has critical ground-state correlations, and 
correspondingly gapless collective excitations [51 ITTl IT2"]. 
However, when the dynamics are generalized to generate 
even a small concentration of second neighbor dimers, 
the ground-state correlations become short-ranged and a 
gap opens in the spectrum. Indeed, although the ground- 
state correlations in this limit can no longer be computed 
directly using Pfaffian methods|13j. close to the critical 
point, where the concentration of second neighbor dimers 
is small, the correlations are asymptotically equivalent to 
those of a 1+1 dimensional massive Thirring model, and 
so are still exactly known [Tl]. 



In the course of this study, we have also recognized 
a new, possibly useful diagnostic tool for numerical 
searches for spin-liquids: We find the exact spin-liquid 
ground state on a torus with a finite circumference, L y , 
spontaneously breaks translational symmetry (forms a 
"columnar" density wave state) for odd L y , but is trans- 
lationally invariant for even L y . (A similar phenomenon 
occurs when a fractional quantum Hall fluid is considered 
in the narrow torus limit |15U16j .) For odd L y , the ampli- 
tude of the columnar order parameter decays exponen- 
tially in proportion to —jjj ex P[~ L y /2£\, where £ is the 

dimer-dimer correlation length[17]. We note that a simi- 
lar even-odd effect with exponential decrease in the am- 
plitude of the columnar order parameter for odd L y have 
been observed[3] in DMRG studies of the square-lattice 
spin-1/2 J\-Ji model on cylinders with L y = 3 — 10. 

Model: The quantum dimer model is defined on a 
Hilbert space with distinct, orthonormal states corre- 
sponding to each allowed hard-core dimer covering of 
the lattice. The Hamiltonian is defined by matrix el- 
ements between dimer states, with "potential" terms, 
V, V, XV", and V" /A, which are diagonal in the dimer 
basis and associate an interaction energy with various lo- 
cal arrangements of dimers, and "kinetic" terms, t, if , 
and t" , which involve a local rearrangement of a small 
number of dimers. On the square lattice, we represent 
the Hamiltonian graphically as: 
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where black bonds are occupied by dimers, the first 
sum runs over all plaquettes, and the second and third 
over all pairs of adjacent plaquettes (with both orienta- 
tions). The first line is Eq. |l]) is the original quantum 
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dimer model on the square lattice with only first neigh- 
bor dimers [9] ■ The added terms are the shortest-range 
terms involving second-neighbor dimers. The parameter 
A determines the relative preference for first and second- 
neighbor dimers. 

A similar construction can be used to define the model 
on the honeycomb lattice: 
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where various parameters have similar meaning as in the 
case of square lattice, and it is implicit that all terms are 
to be summed over symmetry related orientations. 

In defining these models, we have restricted the Hilbert 
space by excluding configurations with crossed dimers. 
We do this for two reasons: Firstly, more microscopi- 
cally, the dimers are thought to be representations of spin 
singlets (or valence bonds) which are assumed to be the 
building blocks of the low energy subspace of an under- 
lying quantum spin 1/2 problem with strong frustration. 
However, there are only two linearly independent sin- 
glet states corresponding to four spin 1/2 's on a plaque- 
tte, or in other words, if we identify dimer configurations 
with orthogonalized versions of valence bond states, then 
|X) oc — (|I I) + |")). (It is an open question whether 
the remaining states with first and non-crossing sec- 
ond neighbor valence bonds are linearly indep endent [ 1 8] . ) 
Secondly, for the model as defined, crossed dimers would 
be non-dynamical unless we were to include additional 
terms in the Hamiltonian. 

For the system with a torus or cylinder geometry, it 
is known that dimer configurations can be classified ac- 
cording to topological labels or winding numbers; dimer 
configurations with different winding numbers cannot be 
connected by local moves of dimers. On the square 
and honeycomb lattices, configurations with only first 
neighbor dimers are characterized by integer winding 
numbers[H] due to the bipartiteness of the model. How- 
ever, for the dimer configurations allowing both first and 
second neighbor dimers, the bipartiteness is lost and on 
a 2-torus there are only four topological sectors labeled 
by W = (W x ,W y ), where W, is odd or even. A simple 
way to determine W x is to draw a vertical line across the 
whole system and see how many dimers it cuts whose 
parity is then W x , and W y is defined similarly in terms 
of the even or oddness of the number of dimers crossing 
a horizontal line. More generally, the number and char- 
acter of the topological sectors depends on the boundary 
conditions in ways that are straightforward to determine. 

The RK point: For generic parameters the exact 
ground state is not known. However, at the so-called 
Rohksar-Kivelson (RK) points, exact ground states of 



both the models can be determined explicitly. For both 
Hamiltonians [Eq. ([!]) and Eq. p])], the generalized RK 
points are given by 



t = V, t' = V'., and t" = V", 



(3) 



and any finite A. At these RK points, both Eq. ([I]) and 
Eq. ^ can be expressed as a sum of projection operators 
and can thus be shown to be positive semi-definite [9], 
with groundstates: 
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where the summation is over all dimer configurations, 
c, in a given topological sector and n s (c) is the number 
of second-neighbor dimers in c. That this is an exact 
ground-state can be seen by explicitly checking that it is 
anihilated by the Hamiltonian [Eq. ([I]) or Eq. Q]. On 
the torus, there are four sectors so the ground states have 
fourfold topological degeneracy [TJJ [2D]. 

It remains to characterize the phases described by 
these exact ground-states in terms of the behavior of 
the dimer-dimer correlation functions. For the RK wave 
function with only first neighbor dimers, the ground state 
is known to correspond to a quantum multi-crtical point 
with power-law decay of correlation functions. However, 
as we shall show, the wave function described by Eq. Q 
is more akin to that of the dimer model on the triangu- 
lar lattice in the sense that all groundstate correlation 
functions fall exponentially with distance[21 . 

Path integral representation: The dimer-dimer 
correlation function in a ground state RK wave function 
| "Jo) is defined as 
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where Dij denote dimer operators on the link (ij) which 
is 1 when the link is occupied by a dimer and oth- 
erwise. This is equivalent to the correlations of a clas- 
sical dimer model in which the wave function normal- 
ization, Z = (v^ol Wq), plays the role of the partition 




FIG. 1. The arrow pattern on on the square and honeycomb 
lattices used in defining the Grassmann path integral repre- 
sentation of the classical dimer model. The black and white 
sites are the two sublattices. 
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function. For wave functions involving only first neigh- 
bor dimers, it is known that dimer densities and dimer- 
dimer correlation functions can be computed in terms of 
a path integral representation of non-interacting Grass- 
mann fields, since a consistent arrow pattern determining 
the action of Grassmann fields is always possible for a pla- 
nar graph [T^]. Due to the appearance of second neighbor 
dimers, the graph is no longer planar and the Pfaffian 
method does not work directly. Nonetheless, the parti- 
tion function can still be exactly expressed by the follow- 
ing path integral representation of interacting Grassmann 
fields (see the proof in Appendix) : 
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where are Grassmann numbers and [da] — Jli!=i den 
(N is the number of sites), tij = ±1 on the first neigh- 
bor links, Uj = ±A on second neighbor links, and (ijkl) 
denotes a pair of crossing second neighbor links ij and 
kl. Here, the sign of tij is determined from the pattern 
of arrows shown in Fig. [Ij such that Uj is positive if the 
arrow on the link (ij) points from i to j, and negative if 
it points from j to i. (The choice of a particular pattern 
of arrows amounts to a choice of gauge.) 

Since the new thing here is the four-field "interaction" 
term, it worth discussing its origin heuristically. In the 
absence of the interaction term, the path integral has 
a contribution from dimer configurations with crossed 
dimers with an associated amplitude factor, [— t^t^i], 
which can be negative. The interaction term likewise 
generates crossed dimers with an amplitude tijtki, which 
is just what is needed to cancel the spurious term in the 
absence of interactions. The proof is in Appendix. 

Massive Thirring model: The path integral in 
Eq. ^ describes the classical dynamics of fermionic fields 
in 2 spatial dimensions which can be mapped to a 1+1 
dimensional quantum theory by treating one spatial di- 
mension, say y, as time t. 

For the square lattice, the arrow pattern we have cho- 
sen doubles the unit cell; for the honeycomb lattice, the 




FIG. 2. The Brillouin zone (BZ) for the square and honey- 
comb lattices: The solid lines mark the boundaries of the full 
BZ corresponding to the arrow patterns in Fig. 1, with two 
sites per unit cell. The shaded area is the "half BZ" over 
which the fc-sum in Eq. 9 is carried out. The dot indicates 
the location of the Dirac node, K, in the limit A — > 0. 



unit cell already contains two sites, and no further in- 
crease is necessary in order to define a consistent pattern 
of arrows; we will thus define a new two-component field 
with a pseudo-spin index, a, which labels the sublattice. 
Translation symmetry allows us to block diagonalize the 
non-interacting part of the Hamiltonian by Fourier trans- 
form in terms of a new set of complex Grassmann fields 
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since a 7 - are real Grassmann fields. 



The non-interacting part of the action is then 
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where the sum runs over only half the Brillouin zone (See 
Fig. 2.) (since we have combined real fields with k and 
—k into a single complex field), and 
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for the square and honeycomb lattice, respectively. Here 
ij are the vectors connecting first neighbors of the hon- 
eycomb lattice: e.\ = (0,-1), e 2 = (v^/^, 1/2), e 3 = 
(-■v/3/2, 1/2), and e j+3 = ej. 

When A = 0, the spectrum of is gapless and has 
a single Dirac point at a momentum labeled K (the 
dot points in Fig. [2]) within the half Brillouin zone: 
K = (0,7r/2) for the square and K = (4tt/3\/3,0) for 
the honeycomb lattice. A finite A opens up a gap in the 
spectrum and is the mass term of the Dirac fermions. 
Indeed, taking the continuum limit of Eq. [9] for small A 
produces a formal equivalence to the 1+1 D Dirac action: 

Sq=vf / dxdy^ [id y (T x — id x a y + m cr 2 ]V>, (11) 



where for the square lattice, vf — 2 and rrio = 2 A while 
for the honeycomb lattice, vp = I and mo = 2\/3A. 
Here ipa(r) is the slowly varying piece of the inverse 
Fourier transform of b^ , with the rapidly oscillating fac- 
tor, exp[iK ■ r\, removed. For small A, the interaction 
term can also be treated in the naive continuum limit, so 



(12) 
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where g = 8A 2 for the square lattice and g — 24A 2 for the 
honeycomb is the strength of attractive interaction. The 
effective interaction is attractive as it serves to cancel the 
(negative) weight for crossed dimers generated by Sq. 
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Rescaling the spatial coordinates so that vf = 1, iden- 
tifying y with the imaginary time, t, and introducing the 
notation ip = ip^a z , we see that Eq. (12) is precisely the 
Euclidean action of the massive Thirring model. The 
massive Thirring model is, in turn, exactly solvable by 
Bethe ansatz[T3]. The mass is renormalized as: 

m^m e- 9/vF7r . (13) 

In other words, the mass gap gets renormalized down but 
survives under weak attractive interactions. Indeed, the 
effect of interactions is only perturbative and does not 
qualitatively change the feature of the non-interacting 
model, as long as A is small, i.e., A 2 log(-^) <C 1. 

The finite gap in the spectrum at small but finite A im- 
plies that the ground state dimer correlations are short- 
ranged [22 j. Consequently, the ground states are gapped 
Z2 quantum spin liquids with a fourfold topological de- 
generacy on a torus [23] • This is our central result. 

Even-odd effect: Unambiguously identifying gapped 
spin-liquid phases is notoriously difficult experimentally 
since spin-liquid ground states exhibit no conventional 
broken symmetries. In this context, we note that when 
a system with a gapped spin-liquid groundstate is placed 
on a torus with a finite odd-length circumference L y in 
one direction, a two-fold degenerate density-wave state 
(i.e. which spontaneously breaks translational symme- 
try) necessarily arises. However, the density wave or- 
der vanishes when L y is even. This even-odd effect 
is rooted in the seminal theorem of Lieb, Schultz and 
Mattis[24l [25] : For a ID system with an odd number 
of electrons per unit cell, there must always exist a dis- 
tinct state with momentum tt relative to the ground-state 
whose energy approaches arbitrarily close the ground 
state energy in the thermodynamic limit. (This theorem 
can be applied to the quantum dimer model by associat- 
ing each dimer with a singlet pair of electrons [9].) Ap- 
plying the Lieb-Schultz-Mattis theorem to gapped spin- 
liquids, one concludes that this ground-state degeneracy 
must correspond to translational symmetry breaking for 
L y odd. 

We have calculated the density-wave order parameters 
for various quantum dimer models on a torus of finite 
circumference L y . In all cases, the density- wave order 
vanishes for even L y . However, for odd and large L y , we 

find columnar order which decays as L y exp[— L y /2£], 
where £ is the dimer-dimer correlation length in 2D. For 
the square lattice with small A, the columnar density- 
wave order parameter C(L y ) (defined as the difference 
between strong and weak bonds) is 



C(L y ) 



/4A(1 -4A 2 ) 
tt(1-8A 2 ) 



-i)/2? 
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for L y 3> £, where £ ~ 1/(4A). Both the even-odd effect 
and exponential decay of the columnar order with in- 
creasing L y have been observed in DMRG studies of the 
spin-1/2 J\-Ji Heisenberg model on a square lattice[3J. 



We have obtained similar even-odd results for the mod- 
els on the triangular and Kagome lattices. Because the 
honeycomb lattice has two sites per unit cell, it is nec- 
essary to introduce twisted boundary conditions in order 
to have a torus with odd L y . This twist globally lifts the 
degeneracy between the two, distinct columnar states. 
None-the-less, the magnitude of the bond alternation de- 
cays exponentially with increasing L y in this case, as well. 

Remark: The present analysis of the ground-state 
correlations applies only for A small; it is uncertain 
whether or not there is a transition, possibly to a broken 
symmetry state, above a nonzero critical A. 
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SUPPLEMENTAL MATERIAL: PROOF OF THE 
PATH INTEGRAL REPRESENTATION 

Here, we prove that the path integral over interacting 
Grassmann fields produces an exact representation of the 
Boltzman weight for the classical dimer problem defined 
by the ground-state expectation- values at the generalized 
RK points on the square and honeycomb lattices. For 
simplicity, to avoiding the complication arising from the 
existence of different topological sectors for the system on 
a torus, we will first treat the system with open boundary 
conditions along both direction. We turn to the problem 
on a torus subsequently. Moreover, we first focus on the 
calculation of the partition function of the classical dimer 
problem (the wave-function normalization), and turn to 
the calculation of ground-state correlation functions in 
the final paragraphs. 

Define zq to be the weighted sum of legitimate dimer 
configurations without crossed dimers which determines 
the self-overlap of the ground state RK wave function 
l*o>: 



z = (*o|*o>=£A"< 
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We further define a set of "illegitimate" partition func- 
tions, z n , by the same expression, but with the sum over 
configurations corresponding to all dimer coverings with 
exactly n > pairs of crossed dimers. The partition 



function from the path integral is Z 
goal is to prove that Z = z . 
It is directly clear that 

N d /2 

z = Z m , 
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where Z m is the m th term in the expansion of Z in powers 
of the interaction strength (i.e. the second term in Eq. 
7), and N d is the number of dimers or one half the number 
of lattice sites N d = N/2 (which we assume to be even 
for simplicity). Higher order terms, Z m with m > Nd/2, 
vanish identically due to the nature of the Grassmann 



integral. For Z , we obtain 
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where the sign (—1)™ is due to the "incorrect" arrow 
pattern involving crossed pairs of dimers generated by 
the non-interacting term in S, and the factor l/N d \ is 
canceled by the expansion of (—S) Nd in terms of dimer 
configurations. 

The first term involving interactions is 



Zx = \ [da] 
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= zi - 2z 2 H (-1) 2 I 1 jz Nd/2 
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71=1 ^ ' 

where the coefficients (™) above is due to the fact that 
there are (™) number of ways to choose a crossed pla- 
quette generated by the interaction term. In general, we 
obtain 

/(_ a\N d -m 
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By denoting Z 



1lZ=o fnZ n , we finally obtain 
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In other words, we have proved that Z = z . 

The path integral representation above is exact for a 
system with open boundary conditions along both direc- 
tions. For a system on torus, the dimer configurations 
break into topological sectors, W — (W x ,W y ) where 
W a can take on values e or o depending on whether 
the number of dimers crossing a reference line in the 
direction perpendicular to e a is even or odd. The RK 
wave function in a given topological sector is given by a 
combination of path integrals of Grassmann fields with 
different boundary condition. Let z\y denote the RK 
wave function self-overlap in sector W. For W — (e, e), 
Z(e,e) = {Z pp + Z ap + Z pa + Z aa ) I '4, where p/a means that 
the path integral is computed for Grassmann fields with 
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pcriodic/anti-periodic boundary conditions along x and y 
directions. Similarly, z^ a ) = (-Z pp + Z pa - Z ap + Z aa ) / '4, 

z {o,e) = {~Zpp — "^pa + ^ap + Z aa )/A, and Z(o,o) — 
(~Zpp + ^pa + ^ap — Zaa)/ 4- 

To compute dimer densities and dimcr-dimer corre- 
lation functions in the RK wave function, we use the 
same path-integral representation of interacting Grass- 
mann fields. The dimer density on a given link (ij) is 

(Aj) ^*o|A;l*o> 
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where sgn is the sign function which is used to make 
the density non-negative. The explicit expressions are 
slightly different between first neighbor and second neigh- 
bor links. For first neighbor links, the density density is 
given by 



(A,) = 



sgn(ijj) 
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The density density on second neighbor links (ij) is 



sgn(ty) 



[da] (-iaidj)(l- ^ itkiakcn) 
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Similarly, we can obtain the dimcr-dimer correlation 
functions from the path integral representation. The 
dimcr-dimer correlation function is given by 



(*o| *o> 
1 , . 5 2 Z 
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for which an explicit expression in terms of the Grass- 
mann path integral can be written out straightforwardly. 



